/* Geometry.h
 *
 *	Copyright (C) 2008 CRIMERE
 *	Copyright (C) 2008 Jean-Marc Mercier
 *	
 *	This file is part of OTS (optimally transported schemes), an open-source library
 *	dedicated to scientific computing. http://code.google.com/p/optimally-transported-schemes/
 *
 *	CRIMERE makes no representations about the suitability of this
 *	software for any purpose. It is provided "as is" without express or
 *	implied warranty.
 *
 *  This program is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *   You should have received a copy of the GNU General Public License
 *   along with this program.  If not, see <http://www.gnu.org/licenses/>.
 *
 *
 */
#if !defined(_Operators_)
#define _Operators_
/*! \file Operators.h
    \brief Arborescent classification of Operators on data.

	This file aims to classify operators acting on data set to provide basic adapted algebra operations
	on operators.
*/
#include <boost/function.hpp>
#include <src/miscellaneous/io/IOFacilities.h>
#include <src/miscellaneous/OS_basic.h>
#include <src/data/Matrix.h>
/*! \brief helper for finding outer multiplication and division type.

	the outer type for multiplication and division is called  base_value_type.
	- if A is of fundamental type (int, float, long, etc..), then we suppose that the external product coincide with the inner product and base_value_type is A.
	- if A is not of fundamental type, then A must define its outer product and declare its outer type in A::base_value_type.
*/
template <class A,class R>
struct outer_structure
{
	typedef boost::integral_constant<bool, 
      ::boost::is_fundamental<A>::value > truth_type;

	template <class T>
	struct outer_helper{
		typedef typename A::base_value_type type;
	};
	template <>
	struct outer_helper<OS_true>{
		typedef typename A type;
	};
//! base_value_type
	typedef typename outer_helper<truth_type>::type base_value_type; // default is A for integral type, A::value_type for the others
};


/*! \brief	common naming base class for finite difference operators
*/
struct Finite_Diff_Operator {
//!return " Fin. Diff. "
	static std::string get_name() { return " Fin. Diff. ";};
};
/*! \brief	common naming convention for particle<OS_double> operators
*/
struct Particle_Operator {
//!return " Particle "
	static std::string get_name() { return " Particle ";};
};

/*! \brief interface class for functions.
	This class is the interface class for unary function.
	it shows an	result_type operator()(argument_type) const that trow an excpation.
*/
template<class _Arg, class _Res>
struct operator_base : public std::unary_function<const _Arg& , _Res >
{
	typedef typename _Arg				argument_base_type; 
	typedef operator_base				type;
	typedef type						Base_Cast;
	typedef smart_ptr<type>				smart_ptr_function;
//!	Throw an exception at run time : must be overloaded in derived classes.
	virtual result_type operator()(argument_type data) const {
		OS_DYNAMIC_CHECK(false,"operator_base : this virtual must be overloaded in derived class");
		return result_type();
	}; 
};
/*! \brief interface class for binary functions.
	This class is the interface class for binary function.
	owns a binary operator that trow an exception at run time.
*/
template<class _Arg, class _Res>
struct operator_binary_base : public std::binary_function<const _Arg& ,const _Arg&,  _Res >
{
	typedef typename _Arg argument_base_type; 
//!	Throw an exception at run time : must be overloaded in derived classes.
	virtual result_type operator()(first_argument_type left,second_argument_type right) const {
		OS_DYNAMIC_CHECK(false,"operator_binary_base : this virtual must be overloaded in derived class");
		return result_type();
	};
};
/*! \brief	derivable function interface
*/
template <class data, class result = data>
struct derivable : public virtual operator_base<data,result>
{
	typedef operator_base<data,result>			Base;
//!	Throw an exception at run time : must be overloaded in derived classes.
	virtual smart_ptr_function get_derivative() const {
		OS_DYNAMIC_CHECK(false,"integrable : this virtual must be overloaded in derived class");
		return smart_ptr_function();
	};
};
/*! \brief	integrable function interface
*/
template <class data, class result = data>
struct integrable : public virtual operator_base<data,result>
{
	typedef operator_base<data,result>			Base;
//!	Throw an exception at run time : must be overloaded in derived classes.
	virtual smart_ptr_function get_integral() const {
		OS_DYNAMIC_CHECK(false,"integrable : this virtual must be overloaded in derived class");
		return smart_ptr_function();
	};
};
/*!
	This class defines the behavior 2) of the class Operator.
	It defines a graph of connectivity between data set.
	Due to the choice of the graph, we provide the cardinal function.
	maybe a choice ok kind std::map<indice_type, std::set<indice_type>> would be
	more pertinent.
*/
template<class A
			, class Data = weight_set_default_def<typename A::base_value_type> 
			, class graph_type = std::set< std::pair<typename A::indice_type, typename A::indice_type> > >
struct graph_base
{
	typedef Data																	weight_set_type;
	typedef typename A::indice_type													indice_type;
	typedef typename graph_type::value_type											edge_type;
	typedef std::size_t																size_type;
	graph_base(const weight_set_type & from) : weight_set_(from)
	{};
	graph_base() {};
	virtual size_type cardinal(const indice_type & i) const {
		OS_DYNAMIC_CHECK(false,"to implement for your own data set");
		return 0;
	};
	virtual const graph_type& get_graph() const{
		OS_DYNAMIC_CHECK(false,"to implement for your own data set");
		return graph_type();
	};
	weight_set_type weight_set_;
};
/*!Specialization of graph_base for the default template paramaters.

*/
template<class A>
struct graph_base<A
	, weight_set_default_def<typename A::base_value_type>
	, std::set< std::pair<typename A::indice_type, typename A::indice_type> > >
{
	typedef weight_set_default_def<typename A::base_value_type>				
																					weight_set_type;
	typedef std::set< std::pair<typename A::indice_type, typename A::indice_type> >
																					graph_type;
	typedef typename A::indice_type													indice_type;
	typedef typename graph_type::value_type											edge_type;
	typedef std::size_t																size_type;
	graph_base(const weight_set_type & from) : weight_set_(from){};
	graph_base() {};
	virtual size_type cardinal(const indice_type & i) const {
		OS_DYNAMIC_CHECK(false,"to implement ");
		return 0;
	};
	virtual const graph_type& get_graph() const {
		return graph_;
	};
	virtual void build_graph() {
		typedef weight_set_type::const_iterator1 const_iterator1;
		typedef weight_set_type::const_iterator2 const_iterator2;
		const_iterator1 itl;
		const_iterator2 itl_;
		size_type size = weight_set_.nnz();
		if (graph_.size() != size) {
			graph_.clear();
			const_iterator1 end(weight_set_.end1());
			for (itl = weight_set_.begin1();itl != end;++itl) {
					indice_type i = itl.index1();
					const_iterator2 end2(itl.end());
					for (itl_ = itl.begin();itl_ != end2;++itl_) {
						indice_type j = itl_.index2();
						graph_.insert(std::make_pair(i,j));
					};
			};
		};
	};

	weight_set_type weight_set_;
	mutable graph_type graph_;
};



/*! state variable pattern.

	This class defines the behavior 3) of the class Operator : each time that a variable of state changes
	(a call to set_state_variable), then the function update() must check that the connectivity graph
	\{(i,j)}_{i \in \Omega_S, j \in \nu_i} is correctly defined.
*/
template <class T>
struct state_variable_traits
{
//!	definition of the state variable type
	typedef typename T				state_variable_type;
	virtual ~state_variable_traits(){};
	state_variable_traits(){};
	const state_variable_type & get_state_variable() const {return state_variable_;};
//!	set the state variable type and call update()
	void set_state_variable(const state_variable_type & from) {
		state_variable_ = from;
		update();
	};

private :
//!	virtual update() function : throw exception at run time. Must be overloaded in derived class.
	virtual void update(){
		OS_DYNAMIC_CHECK(false,"this virtual must be overloaded in derived class");
	}; 
//!	The state variable variable
	state_variable_type state_variable_;
};

/*! \bried Sized decision for update
	
	A common feature for operators : the graph is recomputed only if the size change.
*/
template <class T>
struct sized_variable_traits : public state_variable_traits<T>
{
public :
	typedef  typename state_variable_traits<T>									Base;
	sized_variable_traits() : size_(0) {};
//! pure virtual rebuild function
	virtual void rebuild() = 0;
//! update call rebuild only if the state_variable_size changes. see state_variable.
	void update() {
		if (get_state_variable()->size() != size_) {
			size_ = get_state_variable()->size();
			if (size_ > 0)
				rebuild();
		};
	};
protected :
	OS_size size_;
};

/*! \brief common boundary definition traits
*/
template<class _Arg>
struct boundary_base {
//! Dim is the dimension of the unrgument
	enum :std::size_t{Dim = _Arg::Dim};
//! Boundary type tag structure (eg Periodic, Dirichlet,...)
	typedef typename _Arg::Boundary Boundary;
};


/*! \brief Base class for all operator. This class cannot be instanciated.

	This class defines the semantic of Operators.
	Let f be an instance of Operator.
	It must have the following behavior :
	-# It define a function of argument S_type and return an object of kind R_type.
	It means that we can write
	R_type R = f(S);
	basically, this behavior is decribed by the structure operator_base. Thus the semantic
	Operator::argument_type, Operator::result_type is provided.
	-# The input data S describes a collection of data \f$S = \{S_i\}_{i \in \Omega_S}\f$, 
	and return a collection of data \f$R = \{R_i\}_{i \in \Omega_R}\f$, where 
	\f$\Omega_S\f$, \f$\Omega_R\f$, are finite set of indices.
	To that aim, they must define the connectivity between each data, that is they acts as
	\f$R_i = f( \{S_j\}_{j \in \nu_i}\f$, where \f$\nu_i\f$ is a subset of 
	\f$\Omega_S\f$. In other words, each data \f$S_i\f$ is connected
	to the data \f$\{S_j\}_{j \in \nu_i}\f$. This mechanism defines a graph
	 \f[
		\ {i,j}_{i \in \Omega_S, j \in \nu_i}
	\f]
	,and an operator must provide this graph structure.
	In term of container, we chose to define this graph as a structure of kind std::set<std::pair<indice_type, indice_type>>, where indice_type is deduce 
	within the class S (see the structure graph_base).
	-# Operator depends on an external variable to compute themselves. This extra variable is called 
	a variable state.
	The basic mechanism is described by the structure state_variable_traits :
	Each time that this state variable is updated, a call to the function update() to define the graph is achieved.
	Note : we tried to embed this behavior within an Observable / Observer pattern, but it was not adapted.
	-# A boundary tag given by boundary_base. 

	To resume, an operator must provide the following implementation :
	result_type operator()(argument_type) const =0; implementation of the function itself.
*/
template<class _Arg, class _Res, class State_Variable = smart_ptr<_Arg> >
struct Operator : public virtual operator_base< _Arg , _Res >,
	public state_variable_traits<State_Variable>,
	public virtual graph_base<_Arg>,
	public boundary_base<_Arg>
{

	Operator(){};
	virtual ~Operator(){};
};
/* /brief non compilable base class branching for binary operators

	We call a binary operator an Operator that act as \f$ R = f(S) \f$, where
	\f$ R=\{R_i\}_i \f$ and \f$ S=\{S_i\}_i \f$ are collections, and f acts as
	\f[
		R_i = \sum_{j \in \nu_i} \omega_{i,j}f(S_i,S_j)
	\f]
	The connectivity graph depends only on the size of the state variable.
*/
template<class _Arg, class _Res, class State_Variable, class BinaryOperator >
struct BinaryOperator : public Operator<_Arg, _Res,  State_Variable>
{
	OS_STATIC_CHECK(false);
};
/*!
	This class provide a bind for non time dependant Operators.
*/
template<class Op, typename time_type = OS_double>
class BindTimeOperator : 
	public Op
{
public :
	typedef typename Op::argument_type argument_type;
	typedef typename Op::result_type result_type;
	typedef typename Op::argument_base_type argument_base_type;
	typedef typename time_type					time_type;
	BindTimeOperator(const Op & from) : Op_(from) {};
	BindTimeOperator(const BindTimeOperator & from) : Op(from) {};
	BindTimeOperator(){};
	virtual ~BindTimeOperator(){};
	result_type operator()(const time_type & time,typename argument_type arg) const {
		return Op::operator()(arg);
	};
};

template <class Op>
	OS_double inf_distance(const Op & op)
	{
		typedef typename Op												operator_type;
		typedef typename Op::indice_type								indice_type;
		typedef typename operator_type::state_variable_type				state_variable_type;
		typedef typename operator_type::const_iterator					const_iterator;
			OS_double return_(OS_MAX_REAL);
			const state_variable_type & part  = op.get_state_variable();

			if (part.size() > 0) {
				const_iterator = op.begin(), end = op.end();
				for (;beg != end;++beg) {
					indice_type i = beg->first, j = beg->second;
					return_ = std::min(return_, jmmath::norm2(part[i]-part[j]) ) ;
				}
			}
		return return_;
	};







#endif
